// $Id: import-pass2.c 10 2013-09-13 05:40:46Z alan.watson.f@gmail.com $

// Copyright 2013 Alan M. Watson

#include "jatte.h"

#include <stdio.h>
#include <stdlib.h>

static void
read_pgm(FILE *fp, size_t *nyp, size_t *nxp, unsigned int ***zp)
{
  size_t nx;
  size_t ny;
  unsigned long maxz;
  if (fscanf(fp, "P5\n%zd %zd\n%lu\n", &nx, &ny, &maxz) != 3) {
    jatte_error("the input file is not a binary PGM file.");
  }
  
  unsigned int **z = jatte_checked_alloc_2d_unsigned_int_array(ny, nx);
  
  for (size_t iy = 0; iy < ny; ++iy) {
    for (size_t ix = 0; ix < nx; ++ix) {
      int b0 = getc(fp);
      int b1 = getc(fp);
      if (b0 == EOF || b1 == EOF) {
        jatte_error("the input file is truncated.");
      }
      z[iy][ix] = (unsigned int) b0 * 0x100u + (unsigned int) b1;
    }
  }
  
  *nxp = nx;
  *nyp = ny;
  *zp = z;  
}

static void
write_pgm(FILE *fp, size_t ny, size_t nx, unsigned int **z)
{
  fprintf(fp, "P5\n%zd %zd\n65535\n", nx, ny);
  for (size_t iy = 0; iy < ny; ++iy) {
    for (size_t ix = 0; ix < nx; ++ix) {
      unsigned int u = z[iy][ix];
      unsigned int b0 = u / 0x100u;
      unsigned int b1 = u % 0x100u;
      putc(b0, fp);
      putc(b1, fp);
    }
  }  
}

int
main(int argc, char **argv)
{
  if (argc != 2)
    jatte_error("usage: %s fits_prefix", argv[0]);
  const char *fits_prefix = argv[1];

  // Read the raw pixels from the PGM file produced by dcraw on stdin.
  size_t nx;
  size_t ny;
  unsigned int **z;
  read_pgm(stdin, &ny, &nx, &z);
  
  // Extract the B, G0, G1, and R channels.
  unsigned int **b  = jatte_checked_alloc_2d_unsigned_int_array(ny / 2, nx / 2);
  unsigned int **g0 = jatte_checked_alloc_2d_unsigned_int_array(ny / 2, nx / 2);
  unsigned int **g1 = jatte_checked_alloc_2d_unsigned_int_array(ny / 2, nx / 2);
  unsigned int **r  = jatte_checked_alloc_2d_unsigned_int_array(ny / 2, nx / 2);
  for (size_t iy = 0; iy < ny / 2; ++iy) {
    for (size_t ix = 0; ix < nx / 2; ++ix) {
      b [iy][ix] = z[iy * 2 + 0][ix * 2 + 0];
      g0[iy][ix] = z[iy * 2 + 0][ix * 2 + 1];
      g1[iy][ix] = z[iy * 2 + 1][ix * 2 + 0];
      r [iy][ix] = z[iy * 2 + 1][ix * 2 + 1];
    }
  }
  
  // Write the FITS files.
  {
    char *name = jatte_mprintf("%s-b.fits", fits_prefix);
    FILE *fp = jatte_checked_fopen(name, "w");
    jatte_write_fits_2d_unsigned_int(fp, ny / 2, nx / 2, b);
    fclose(fp);
    free(name);
  }
  {
    char *name = jatte_mprintf("%s-g0.fits", fits_prefix);
    FILE *fp = jatte_checked_fopen(name, "w");
    jatte_write_fits_2d_unsigned_int(fp, ny / 2, nx / 2, g0);
    fclose(fp);
    free(name);
  }
  {
    char *name = jatte_mprintf("%s-g1.fits", fits_prefix);
    FILE *fp = jatte_checked_fopen(name, "w");
    jatte_write_fits_2d_unsigned_int(fp, ny / 2, nx / 2, g1);
    fclose(fp);
    free(name);
  }
  {
    char *name = jatte_mprintf("%s-r.fits", fits_prefix);
    FILE *fp = jatte_checked_fopen(name, "w");
    jatte_write_fits_2d_unsigned_int(fp, ny / 2, nx / 2, r);
    fclose(fp);
    free(name);
  }
    
  // Create a W channel.
  unsigned int **w  = jatte_checked_alloc_2d_unsigned_int_array(ny / 2, nx / 2);
  for (size_t iy = 0; iy < ny / 2; ++iy) {
    for (size_t ix = 0; ix < nx / 2; ++ix) {
      w [iy][ix] = b[iy][ix] + g0[iy][ix] + g1[iy][ix] + r[iy][ix];
    }
  }

  // Write the W channel as a PGM file to stdout for further processing.
  write_pgm(stdout, ny / 2, nx / 2, w);

  jatte_free_2d_unsigned_int_array(ny, nx, z);
  jatte_free_2d_unsigned_int_array(ny / 2, nx / 2, b);
  jatte_free_2d_unsigned_int_array(ny / 2, nx / 2, g0);
  jatte_free_2d_unsigned_int_array(ny / 2, nx / 2, g1);
  jatte_free_2d_unsigned_int_array(ny / 2, nx / 2, r);
  jatte_free_2d_unsigned_int_array(ny / 2, nx / 2, w);
  
  return 0;
}